% Estimating confidence three variable VAR with confidence first
%clear all;
%close all;

%load data2.txt
%load ics.txt
%load pfe.txt
%e5y = data2(:,1);
%lrnds = data2(:,2);
%gdp = data2(:,3);
%ics = ics;
%pfe = pfe;

nconsta = 1; %constant term
bps=0; % bandpass filer
tr=200; % bpass long cycle periods
%load var_monte_data;
%macro=var_monte_data;


%macro=xlsread('macrodata_SW_tfp_hours2.xls',1,'B19:R198');
macro=xlsread('macrodata_SW_tfp_hours2_restat.xls',1,'B7:R198');
%macro=xlsread('macrodata_SW_tfp_hours_med.xls',1,'B19:R198');


%macro=xlsread('macrodata_SW_tfp_hours.xls',1,'B19:Q210');  %see xls sheet for details on data: B15 (B19) for the start of 1959q1 (1960q1)
%macro(:,7)=bpass(macro(:,7),2,200);

  macro(:,4)=macro(:,4)/4; %making it to quaterly rate
  macro(:,9)=macro(:,9)/4;%making it to quaterly rate
   %macro(:,9)=log(macro(:,9)/4+1);%making it to quaterly rate
  macro(:,17)=macro(:,17)/4;
  macro(:,16)=macro(:,16)/100;

  
  if bps==1;
 macro(:,4)=bpass(macro(:,4),2,tr);
  macro(:,9)=bpass(macro(:,9),2,tr);
  macro(:,17)=bpass(macro(:,17),2,tr);
  
 macro(:,6)=bpass(macro(:,6),2,tr);   
 macro(:,7)=bpass(macro(:,7),2,tr);
 macro(:,1)=bpass(macro(:,1),2,tr);
 macro(:,2)=bpass(macro(:,2),2,tr);
 macro(:,3)=bpass(macro(:,3),2,tr);
 macro(:,11)=bpass(macro(:,11),2,tr);
 macro(:,16)=bpass(macro(:,16),2,tr);
 macro(:,12)=bpass(macro(:,12),2,tr);
 macro(:,15)=bpass(macro(:,15),2,tr);
 macro(:,14)=bpass(macro(:,14),2,tr);
  end



  
%X=[macro(:,8) macro(:,1) macro(:,2) macro(:,3) macro(:,11) macro(:,4) macro(:,9) ] ; % var_data=20 full(no stock) + no Conf
%X=[macro(:,8) macro(:,1) macro(:,2) macro(:,11) macro(:,4) macro(:,9) macro(:,14) ] ; %var_data=21 no invest(no stock) + Conf

%X=[macro(:,5) macro(:,1) macro(:,2) macro(:,11) macro(:,4) macro(:,9) ]; %%%var_data=22 no invest(no stock) + No Con 2006 TFP 
X=[macro(:,7) macro(:,1) macro(:,2) macro(:,3) macro(:,11) macro(:,9) ] ; % var_data=25 full(no stock) + no Conf +no Interest 2019 TFP
%X=[macro(:,8) macro(:,1) macro(:,2) macro(:,3) macro(:,11) macro(:,9) ] ; % var_data=25 full(no stock) + no Conf +no Interest Adjusted 2019 TFP


%X=[macro(:,7) macro(:,1) macro(:,2) macro(:,3) macro(:,11) macro(:,4) macro(:,9) macro(:,16) ] ; %full
%X=[macro(:,6) macro(:,1) macro(:,2) macro(:,3) macro(:,11) macro(:,4) macro(:,9)   ] ; %with no stock
%X=[macro(:,7) macro(:,1) macro(:,2) macro(:,11) macro(:,4) macro(:,9) macro(:,16)  ] ; % no inv
%X=[macro(:,7) macro(:,1) macro(:,2) macro(:,11) macro(:,4) macro(:,9)   ] ; %no inv &no stock 

%X=[macro(:,7) macro(:,1) macro(:,2) macro(:,3) macro(:,11) macro(:,4) macro(:,9) macro(:,16) macro(:,15) ] ; %full + Conf


%X=[macro(:,7) macro(:,1) macro(:,2) macro(:,3) macro(:,11) macro(:,9) macro(:,14) ] ;%(tfp, y,c,h, pai, conf)
%X=[macro(:,7) macro(:,1) macro(:,2) macro(:,3)  macro(:,9)  ] ;%tfp, y, c, iv,pai
%X=[macro(:,7) macro(:,1) macro(:,2) macro(:,3) macro(:,11)  macro(:,9)  ] ;%tfp, y, c, iv,h,pai
%X=[macro(:,7) macro(:,1) macro(:,2) macro(:,3) macro(:,11) macro(:,4) macro(:,9) macro(:,10) ] ; %tfp, y, c, iv,h,pai, sp




%X=[macro(:,7) macro(:,1) macro(:,2) macro(:,11) macro(:,4) macro(:,9) macro(:,16) ] ; %full(no stock) + Conf

%X=[macro(:,5) macro(:,1) macro(:,2) macro(:,11) macro(:,4) macro(:,9)   ] ;



%X=[macro(:,6) macro(:,12) macro(:,3) macro(:,11) macro(:,4) macro(:,9) macro(:,16)  ] ;

%X=[macro(:,8) macro(:,1) macro(:,12) macro(:,11) macro(:,4) macro(:,9) macro(:,16) ] ;
%X=[macro(:,6) macro(:,1) macro(:,12) macro(:,11) macro(:,4) macro(:,9)] ;

%X=[macro(:,8) macro(:,1) macro(:,12) macro(:,11)] ;

%X=[macro(:,1) macro(:,2) macro(:,3) macro(:,5) macro(:,6) macro(:,7) macro(:,8)] ;

[T,q] = size(X);
p = 4; % Number of lags

%[A,SIGMA,U,V,beta_ols,beta_ols_,A_or,X_or]=olsvarc_hiro(X,p,nconsta);
%[A,SIGMA,U,V,beta_ols,beta_ols_,A_or,X_or,Y_or]=olsvarc_hiro(X,p,nconsta); % original one

% follwoing three comands are for OLS of Kurman, that does not use a
% companion form of VAR to estimate params.
% Using the companion form makes the estimate unreliable due to a large
% inverse matrix to be estimated.

%%%%%%% Not using a companion form of VAR to avoid a large inverse matrix%%
[y, x, vars, lev]=dataset_quarterly_sw(esty1,estq1,esty2,estq2,nlags,shock_extract,use_yields,slope); 
[T,nvars]=size(y);
[A,SIGMA,U,V,beta_ols] = olsvarc_hiro_ols(X, p,nconsta, x,y, nvars, nlags,T);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


B = chol(SIGMA(1:q,1:q));


B = B';

h = 40;

IRF = IRFVAR(A,B,p,h);
IRF0 = IRF;
%[DECOMP,TOT] = vardecomp_eric(IRF,h,p,q);

R = 3000;
%R = 1000; irrelevant?
RR = 3000;
%RR=1000;
pp = 16.50; %16.5 &0.5

nconsta = 0; 
%IMPORTANT: do not include a constant teerm at the bootstrap OLS stage (but include the first regression) 
%since the news-shock identification, i.e., "tfpews.m" uses demeaned OLS (see "tfpnes.m"). 

[BIAS,EE_,EE_1,BIAS_cnt,zz] = boot_firststage_hiro(A,U,q,p,T,V,h,IRF,B,X,R,nlags,nconsta);
[AA, VV_, beta_boots, AA_old,AA_junk1] = boot_second_hiro(A,U,q,p,T,V,h,IRF,B,BIAS,X,RR,pp,nlags,nconsta,BIAS_cnt);

%[BIAS,BIAS_cnt] = boot_firststage_hiro(A,U,q,p,T,V,h,IRF,B,X,R,nlags,nconsta);
%[AA, VV_, beta_boots, AA_old,AA_junk1] = boot_second_hiro(A,U,q,p,T,V,h,IRF,B,BIAS,X,RR,pp,nlags,nconsta,BIAS_cnt);

%AA_tr=[AA VV_]';
%bnew_=AA_tr(:,1:q);


